Housing Prices
Ozan Aygun
11/26/2017
1 Introduction and summary
2 Loading the data
train <- read.csv("train.csv"); test <- read.csv("test.csv")3 Summarizing the data
summary(train)## Id MSSubClass MSZoning LotFrontage
## Min. : 1.0 Min. : 20.0 C (all): 10 Min. : 21.00
## 1st Qu.: 365.8 1st Qu.: 20.0 FV : 65 1st Qu.: 59.00
## Median : 730.5 Median : 50.0 RH : 16 Median : 69.00
## Mean : 730.5 Mean : 56.9 RL :1151 Mean : 70.05
## 3rd Qu.:1095.2 3rd Qu.: 70.0 RM : 218 3rd Qu.: 80.00
## Max. :1460.0 Max. :190.0 Max. :313.00
## NA's :259
## LotArea Street Alley LotShape LandContour
## Min. : 1300 Grvl: 6 Grvl: 50 IR1:484 Bnk: 63
## 1st Qu.: 7554 Pave:1454 Pave: 41 IR2: 41 HLS: 50
## Median : 9478 NA's:1369 IR3: 10 Low: 36
## Mean : 10517 Reg:925 Lvl:1311
## 3rd Qu.: 11602
## Max. :215245
##
## Utilities LotConfig LandSlope Neighborhood Condition1
## AllPub:1459 Corner : 263 Gtl:1382 NAmes :225 Norm :1260
## NoSeWa: 1 CulDSac: 94 Mod: 65 CollgCr:150 Feedr : 81
## FR2 : 47 Sev: 13 OldTown:113 Artery : 48
## FR3 : 4 Edwards:100 RRAn : 26
## Inside :1052 Somerst: 86 PosN : 19
## Gilbert: 79 RRAe : 11
## (Other):707 (Other): 15
## Condition2 BldgType HouseStyle OverallQual
## Norm :1445 1Fam :1220 1Story :726 Min. : 1.000
## Feedr : 6 2fmCon: 31 2Story :445 1st Qu.: 5.000
## Artery : 2 Duplex: 52 1.5Fin :154 Median : 6.000
## PosN : 2 Twnhs : 43 SLvl : 65 Mean : 6.099
## RRNn : 2 TwnhsE: 114 SFoyer : 37 3rd Qu.: 7.000
## PosA : 1 1.5Unf : 14 Max. :10.000
## (Other): 2 (Other): 19
## OverallCond YearBuilt YearRemodAdd RoofStyle
## Min. :1.000 Min. :1872 Min. :1950 Flat : 13
## 1st Qu.:5.000 1st Qu.:1954 1st Qu.:1967 Gable :1141
## Median :5.000 Median :1973 Median :1994 Gambrel: 11
## Mean :5.575 Mean :1971 Mean :1985 Hip : 286
## 3rd Qu.:6.000 3rd Qu.:2000 3rd Qu.:2004 Mansard: 7
## Max. :9.000 Max. :2010 Max. :2010 Shed : 2
##
## RoofMatl Exterior1st Exterior2nd MasVnrType MasVnrArea
## CompShg:1434 VinylSd:515 VinylSd:504 BrkCmn : 15 Min. : 0.0
## Tar&Grv: 11 HdBoard:222 MetalSd:214 BrkFace:445 1st Qu.: 0.0
## WdShngl: 6 MetalSd:220 HdBoard:207 None :864 Median : 0.0
## WdShake: 5 Wd Sdng:206 Wd Sdng:197 Stone :128 Mean : 103.7
## ClyTile: 1 Plywood:108 Plywood:142 NA's : 8 3rd Qu.: 166.0
## Membran: 1 CemntBd: 61 CmentBd: 60 Max. :1600.0
## (Other): 2 (Other):128 (Other):136 NA's :8
## ExterQual ExterCond Foundation BsmtQual BsmtCond BsmtExposure
## Ex: 52 Ex: 3 BrkTil:146 Ex :121 Fa : 45 Av :221
## Fa: 14 Fa: 28 CBlock:634 Fa : 35 Gd : 65 Gd :134
## Gd:488 Gd: 146 PConc :647 Gd :618 Po : 2 Mn :114
## TA:906 Po: 1 Slab : 24 TA :649 TA :1311 No :953
## TA:1282 Stone : 6 NA's: 37 NA's: 37 NA's: 38
## Wood : 3
##
## BsmtFinType1 BsmtFinSF1 BsmtFinType2 BsmtFinSF2
## ALQ :220 Min. : 0.0 ALQ : 19 Min. : 0.00
## BLQ :148 1st Qu.: 0.0 BLQ : 33 1st Qu.: 0.00
## GLQ :418 Median : 383.5 GLQ : 14 Median : 0.00
## LwQ : 74 Mean : 443.6 LwQ : 46 Mean : 46.55
## Rec :133 3rd Qu.: 712.2 Rec : 54 3rd Qu.: 0.00
## Unf :430 Max. :5644.0 Unf :1256 Max. :1474.00
## NA's: 37 NA's: 38
## BsmtUnfSF TotalBsmtSF Heating HeatingQC CentralAir
## Min. : 0.0 Min. : 0.0 Floor: 1 Ex:741 N: 95
## 1st Qu.: 223.0 1st Qu.: 795.8 GasA :1428 Fa: 49 Y:1365
## Median : 477.5 Median : 991.5 GasW : 18 Gd:241
## Mean : 567.2 Mean :1057.4 Grav : 7 Po: 1
## 3rd Qu.: 808.0 3rd Qu.:1298.2 OthW : 2 TA:428
## Max. :2336.0 Max. :6110.0 Wall : 4
##
## Electrical X1stFlrSF X2ndFlrSF LowQualFinSF
## FuseA: 94 Min. : 334 Min. : 0 Min. : 0.000
## FuseF: 27 1st Qu.: 882 1st Qu.: 0 1st Qu.: 0.000
## FuseP: 3 Median :1087 Median : 0 Median : 0.000
## Mix : 1 Mean :1163 Mean : 347 Mean : 5.845
## SBrkr:1334 3rd Qu.:1391 3rd Qu.: 728 3rd Qu.: 0.000
## NA's : 1 Max. :4692 Max. :2065 Max. :572.000
##
## GrLivArea BsmtFullBath BsmtHalfBath FullBath
## Min. : 334 Min. :0.0000 Min. :0.00000 Min. :0.000
## 1st Qu.:1130 1st Qu.:0.0000 1st Qu.:0.00000 1st Qu.:1.000
## Median :1464 Median :0.0000 Median :0.00000 Median :2.000
## Mean :1515 Mean :0.4253 Mean :0.05753 Mean :1.565
## 3rd Qu.:1777 3rd Qu.:1.0000 3rd Qu.:0.00000 3rd Qu.:2.000
## Max. :5642 Max. :3.0000 Max. :2.00000 Max. :3.000
##
## HalfBath BedroomAbvGr KitchenAbvGr KitchenQual
## Min. :0.0000 Min. :0.000 Min. :0.000 Ex:100
## 1st Qu.:0.0000 1st Qu.:2.000 1st Qu.:1.000 Fa: 39
## Median :0.0000 Median :3.000 Median :1.000 Gd:586
## Mean :0.3829 Mean :2.866 Mean :1.047 TA:735
## 3rd Qu.:1.0000 3rd Qu.:3.000 3rd Qu.:1.000
## Max. :2.0000 Max. :8.000 Max. :3.000
##
## TotRmsAbvGrd Functional Fireplaces FireplaceQu GarageType
## Min. : 2.000 Maj1: 14 Min. :0.000 Ex : 24 2Types : 6
## 1st Qu.: 5.000 Maj2: 5 1st Qu.:0.000 Fa : 33 Attchd :870
## Median : 6.000 Min1: 31 Median :1.000 Gd :380 Basment: 19
## Mean : 6.518 Min2: 34 Mean :0.613 Po : 20 BuiltIn: 88
## 3rd Qu.: 7.000 Mod : 15 3rd Qu.:1.000 TA :313 CarPort: 9
## Max. :14.000 Sev : 1 Max. :3.000 NA's:690 Detchd :387
## Typ :1360 NA's : 81
## GarageYrBlt GarageFinish GarageCars GarageArea GarageQual
## Min. :1900 Fin :352 Min. :0.000 Min. : 0.0 Ex : 3
## 1st Qu.:1961 RFn :422 1st Qu.:1.000 1st Qu.: 334.5 Fa : 48
## Median :1980 Unf :605 Median :2.000 Median : 480.0 Gd : 14
## Mean :1979 NA's: 81 Mean :1.767 Mean : 473.0 Po : 3
## 3rd Qu.:2002 3rd Qu.:2.000 3rd Qu.: 576.0 TA :1311
## Max. :2010 Max. :4.000 Max. :1418.0 NA's: 81
## NA's :81
## GarageCond PavedDrive WoodDeckSF OpenPorchSF EnclosedPorch
## Ex : 2 N: 90 Min. : 0.00 Min. : 0.00 Min. : 0.00
## Fa : 35 P: 30 1st Qu.: 0.00 1st Qu.: 0.00 1st Qu.: 0.00
## Gd : 9 Y:1340 Median : 0.00 Median : 25.00 Median : 0.00
## Po : 7 Mean : 94.24 Mean : 46.66 Mean : 21.95
## TA :1326 3rd Qu.:168.00 3rd Qu.: 68.00 3rd Qu.: 0.00
## NA's: 81 Max. :857.00 Max. :547.00 Max. :552.00
##
## X3SsnPorch ScreenPorch PoolArea PoolQC
## Min. : 0.00 Min. : 0.00 Min. : 0.000 Ex : 2
## 1st Qu.: 0.00 1st Qu.: 0.00 1st Qu.: 0.000 Fa : 2
## Median : 0.00 Median : 0.00 Median : 0.000 Gd : 3
## Mean : 3.41 Mean : 15.06 Mean : 2.759 NA's:1453
## 3rd Qu.: 0.00 3rd Qu.: 0.00 3rd Qu.: 0.000
## Max. :508.00 Max. :480.00 Max. :738.000
##
## Fence MiscFeature MiscVal MoSold
## GdPrv: 59 Gar2: 2 Min. : 0.00 Min. : 1.000
## GdWo : 54 Othr: 2 1st Qu.: 0.00 1st Qu.: 5.000
## MnPrv: 157 Shed: 49 Median : 0.00 Median : 6.000
## MnWw : 11 TenC: 1 Mean : 43.49 Mean : 6.322
## NA's :1179 NA's:1406 3rd Qu.: 0.00 3rd Qu.: 8.000
## Max. :15500.00 Max. :12.000
##
## YrSold SaleType SaleCondition SalePrice
## Min. :2006 WD :1267 Abnorml: 101 Min. : 34900
## 1st Qu.:2007 New : 122 AdjLand: 4 1st Qu.:129975
## Median :2008 COD : 43 Alloca : 12 Median :163000
## Mean :2008 ConLD : 9 Family : 20 Mean :180921
## 3rd Qu.:2009 ConLI : 5 Normal :1198 3rd Qu.:214000
## Max. :2010 ConLw : 5 Partial: 125 Max. :755000
## (Other): 9The outcome variable SalePrice is a complete variable, but there are a number of predictors with many missing values. Median house price is $163,000.
par(mfrow = c(1,2))
boxplot(train$SalePrice,main = "SalePrice")
boxplot(log10(train$SalePrice), main = "Log10(SalePrice)")
3.1 Understanding the impact of missing data
Which features have more than 30% missing data?
na.status <- is.na(train)
na.sum <- apply(na.status,2,sum)
names(na.sum) <- colnames(train)
mostly_missing <- which(na.sum > (0.3 * nrow(train)))
na.sum[mostly_missing]## Alley FireplaceQu PoolQC Fence MiscFeature
## 1369 690 1453 1179 14065 variables have substantial amount of missing data. Do they have impact on the house prices?
3.1.1 Alley
Based on the data description: Type of alley access to property
Grvl Gravel
Pave Paved
NA No alley access
library(ggplot2)
ggplot(data = train, aes(x = Alley, y = log(SalePrice), fill= Alley))+
geom_boxplot()
Alley might be important. Gravel homes have significantly lower price. Convert NAs to “NoAlley” to make more sense.
train$Alley <- as.character(train$Alley)
train$Alley[which(is.na(train$Alley))] <- "NoAlley"
train$Alley <- factor(train$Alley)
#Also transform the test set
test$Alley <- as.character(test$Alley)
test$Alley[which(is.na(test$Alley))] <- "NoAlley"
test$Alley <- factor(test$Alley)3.1.2 FireplaceQu
Based on the data description:
FireplaceQu: Fireplace quality
Ex Excellent - Exceptional Masonry Fireplace
Gd Good - Masonry Fireplace in main level
TA Average - Prefabricated Fireplace in main living area or Masonry Fireplace in basement
Fa Fair - Prefabricated Fireplace in basement
Po Poor - Ben Franklin Stove
NA No Fireplaceggplot(data = train, aes(x = FireplaceQu, y = log10(SalePrice), fill= FireplaceQu))+
geom_boxplot()
Let’s look at through a simple linear model:
par(mfrow=c(1,2))
hist(train$SalePrice)
hist(log10(train$SalePrice)) # Looks more gaussian
FirePlaceFit <- lm(log10(train$SalePrice) ~ FireplaceQu, data = train)
summary(FirePlaceFit)##
## Call:
## lm(formula = log10(train$SalePrice) ~ FireplaceQu, data = train)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.38691 -0.09351 -0.01374 0.09604 0.57966
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.50252 0.03013 182.644 < 2e-16 ***
## FireplaceQuFa -0.28729 0.03959 -7.256 9.82e-13 ***
## FireplaceQuGd -0.17957 0.03106 -5.781 1.08e-08 ***
## FireplaceQuPo -0.40442 0.04469 -9.050 < 2e-16 ***
## FireplaceQuTA -0.21003 0.03126 -6.719 3.58e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1476 on 765 degrees of freedom
## (690 observations deleted due to missingness)
## Multiple R-squared: 0.118, Adjusted R-squared: 0.1134
## F-statistic: 25.59 on 4 and 765 DF, p-value: < 2.2e-16par(mfrow=c(2,2))
plot(FirePlaceFit)
It appears that all levels of this feature also have significantly different mean Sales prices. Again, converting NAs to NoFireplace to be more informative.
train$FireplaceQu <- as.character(train$FireplaceQu)
train$FireplaceQu[which(is.na(train$FireplaceQu))] <- "NoFireplace"
train$FireplaceQu <- factor(train$FireplaceQu)
#Also transform the test set
test$FireplaceQu <- as.character(test$FireplaceQu)
test$FireplaceQu[which(is.na(test$FireplaceQu))] <- "NoFireplace"
test$FireplaceQu <- factor(test$FireplaceQu)3.1.3 PoolQC
Based on the data description:
PoolQC: Pool quality
Ex Excellent
Gd Good
TA Average/Typical
Fa Fair
NA No Pool
library(ggplot2)
ggplot(data = train, aes(x = PoolQC, y = log10(SalePrice), fill = PoolQC))+
geom_boxplot()
Note that actually very few houses have pools. Houses with excellent pool condition have significantly higher sales prices. We will also keep this feature, converting NAs to “NoPool”:
train$PoolQC <- as.character(train$PoolQC)
train$PoolQC[which(is.na(train$PoolQC))] <- "NoPool"
train$PoolQC <- factor(train$PoolQC)
#Also transform the test set
test$PoolQC <- as.character(test$PoolQC)
test$PoolQC[which(is.na(test$PoolQC))] <- "NoPool"
test$PoolQC <- factor(test$PoolQC)3.1.4 Fence
Based on the data description:
Fence: Fence quality
GdPrv Good Privacy
MnPrv Minimum Privacy
GdWo Good Wood
MnWw Minimum Wood/Wire
NA No Fence library(ggplot2)
ggplot(data = train, aes(x = Fence, y = log10(SalePrice), fill = Fence))+
geom_boxplot()
FenceFit <- lm(log10(train$SalePrice) ~ Fence, data = train)
summary(FenceFit)##
## Call:
## lm(formula = log10(train$SalePrice) ~ Fence, data = train)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.57500 -0.06747 -0.00520 0.05611 0.72551
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.23634 0.01816 288.272 < 2e-16 ***
## FenceGdWo -0.11851 0.02628 -4.510 9.58e-06 ***
## FenceMnPrv -0.08969 0.02131 -4.210 3.46e-05 ***
## FenceMnWw -0.11305 0.04582 -2.467 0.0142 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1395 on 277 degrees of freedom
## (1179 observations deleted due to missingness)
## Multiple R-squared: 0.08211, Adjusted R-squared: 0.07217
## F-statistic: 8.26 on 3 and 277 DF, p-value: 2.781e-05Interestingly, this feature might have a complex relationship, if any at all, with the response variable. The effects of each of the level deemed as significant, so we will keep this feature at this point as well.
train$Fence <- as.character(train$Fence)
train$Fence[which(is.na(train$Fence))] <- "NoFence"
train$Fence <- factor(train$Fence)
#Also transform the test set
test$Fence <- as.character(test$Fence)
test$Fence[which(is.na(test$Fence))] <- "NoFence"
test$Fence <- factor(test$Fence)3.1.5 MiscFeature
Based on the data description:
MiscFeature: Miscellaneous feature not covered in other categories
Elev Elevator
Gar2 2nd Garage (if not described in garage section)
Othr Other
Shed Shed (over 100 SF)
TenC Tennis Court
NA None
We will deal with the remaining missing value containing features in the next stages of our analysis.
table(train$MiscFeature)##
## Gar2 Othr Shed TenC
## 2 2 49 1ggplot(data = train, aes(x = MiscFeature, y = log10(SalePrice), fill = MiscFeature))+
geom_boxplot()
Again, very few houses have this features, and not sure if these qualities have a real impact on the house value amongst other measured features. Neverthless, I will maintain this feature as well. Just converting NAs to “NoMiscFeature” to be more clear:
train$MiscFeature <- as.character(train$MiscFeature)
train$MiscFeature[which(is.na(train$MiscFeature))] <- "NoMiscFeature"
train$MiscFeature <- factor(train$MiscFeature)
#Also transform the test set
test$MiscFeature <- as.character(test$MiscFeature)
test$MiscFeature[which(is.na(test$MiscFeature))] <- "NoMiscFeature"
test$MiscFeature <- factor(test$MiscFeature)4 Developing expectations from data
Here I will explore the data further, in order to develope a better understanding about the relationships between features and the home prices.
4.1 Continuous features: watching for strong predictors and potential collinearities
Let’s look at the continuous features in the data set:
cont.var <-NULL
for(i in 1:ncol(train)){
if(class(train[,i]) == "integer" | class(train[,i]) == "numeric"){
cont.var <- c(cont.var,i)
}
}4.1.1 Part1: transforming LotArea to make it more gaussian
pairs(log10(train$SalePrice) ~ ., data = train[,cont.var[1:10]], cex = 0.05, pch = 19, col = "navy")
- We need to remove ID, since this is just an identifier for the houses, has no predictive value as expected.
- OverallQual and OverallCond look strong predictors, but might have some collinearity.
- Yearbuilt and YearRemodAdd also look good predictors but most likely contain redundant information.
- For the LOT and SF features, it could be better to look at after log transformation, since the data looks skewed.
par(mfrow=c(3,2))
hist(train$LotArea, breaks = 50, col = "navy");hist(log10(train$LotArea+1),breaks=50,col = "navy")
hist(train$MasVnrArea, breaks = 50,col = "navy");hist(log10(train$MasVnrArea+1),breaks=50,col = "navy")
hist(train$BsmtFinSF1, breaks = 50,col = "navy");hist(log10(train$BsmtFinSF1+1),breaks=50,col = "navy")
Indeed one of these features (LotArea) look more gaussian when it is log10 transformed.
pairs(log10(train$SalePrice) ~ log10(train$LotArea+1) + log10(train$MasVnrArea+1) + log10(train$BsmtFinSF1+1),cex = 0.1, pch = 19, col = "navy")
When transformed, lotArea seems to do a better job in explaining the variability in the SalePrice, other 2 variables however, have substantial 0 values and it is questionable how much this transformation might help.
Therefore, I will only transform the LotArea at this stage and leave the other features as they are, for feature selection steps.
4.1.2 Part2: Area related features generally explain variation in the SalePrice
pairs(log10(train$SalePrice) ~ ., data = train[,cont.var[11:20]], cex = 0.05, pch = 19, col = "navy")
4.1.3 Part3: Size of the garage and number of rooms are good predictors
pairs(log10(train$SalePrice) ~ ., data = train[,cont.var[22:31]], cex = 0.05, pch = 19, col = "navy")
Arguably, these features are all related to the overall area of the house.
4.1.4 Part4: fairly weak predictors
pairs(log10(train$SalePrice) ~ ., data = train[,cont.var[32:38]], cex = 0.05, pch = 19, col = "navy")
Processing the identified features:
# Remove ID variable
train.ID <- train$Id #Hold it in a new variable
train <- train[, -1]
library(dplyr)
train <- dplyr::mutate(train, log10.LotArea = log10(LotArea))
train <- dplyr::select(train, -LotArea)
# Also transform in the test set
test <- dplyr::mutate(test, log10.LotArea = log10(LotArea))
test <- dplyr::select(test, -LotArea)4.2 Categorical predictors: feature engineering to dummy variables
Now let’s look at the categorical variables to see if certain ones are more powerful than others:
cat.var <-NULL
for(i in 1:ncol(train)){
if(class(train[,i]) == "factor"){
cat.var <- c(cat.var,i)
}
}
length(cat.var)## [1] 43We have 43 categorical variables, some of them have only few levels. This is a challenge in many ways. It is not feasible to look each of them individually, therefore, we will first convert them to binary variables and try to use regularization or decision tree-based approaches to select and use most useful features.
library(caret)
# Write a function that gets a data set and converts all factor variables to dummy variables
convert.to.dummy <- function(data.set){
cat.var <-NULL
temp.data <- data.frame(1:nrow(data.set))
for(i in 1:ncol(data.set)){
if(class(data.set[,i]) == "factor"){
cat.var <- c(cat.var,i)
factor.levels <- levels(data.set[,i])
for(j in seq_along(factor.levels)){
dummy.vector = ifelse(data.set[,i] == factor.levels[j],1,0)
dummy.vector <- data.frame(dummy.vector)
colnames(dummy.vector)[1] = paste(names((data.set)[i]),
factor.levels[j],sep = ".")
temp.data <- cbind(temp.data,dummy.vector)
}
}
}
#Remove the original categorical variables from data.set
data.set <- data.set[,-cat.var]
#Add the dummy.variable set
temp.data <- temp.data[,-1] # remove the unnecessary column
data.set <- cbind(data.set,temp.data)
return(data.set)
}
# Process the training set
training.processed <- convert.to.dummy(train)
# Process the test set
test.processed <- convert.to.dummy(test)
# Note that not all the levels are present in the test set, therefore we need to make the features same.
consensus.fearures <- which(colnames(training.processed) %in% colnames(test.processed))
training.processed.SalePrice <- training.processed$SalePrice
training.processed <- training.processed[,consensus.fearures]
training.processed$SalePrice <- training.processed.SalePrice